/*-------------------------------------------------------------------------*/
/**
  @file		contract.c
  @author	Renaud Wanchoor
  @date	July 2004
  @version	$Revision: 2.7 $
  @brief	contract a graph using a computed matching
  
*/
/*--------------------------------------------------------------------------*/
#include "pagrid.h"
/**
 * \brief contract a graph according to a given matching
 * 
 * contract a graph according to the matching graph->match, graph->map
 * too much memory is allocated for the contracted graph and a memory
 * reallocation should be done
 * \param graph graph structure with matching information
 * \param guide permutation used to visit vertices
 * \param cnvtxs number of vertices in coarser graph
 * \return contracted graph structure
 *
 */
GraphType* Contract(GraphType *graph, int *guide, int numVtxs)
{
  int i, j, v, u, k, maxdegree=0;

  GraphType *child;	/*child graph*/
  int *mapadj;		/*storing the label (in the coarser graph)
			  to which each label in the adjncy list 
			  (in the father graph) is mapped*/
  int *record;		/*auxiliary list when searching for adjacent 
			  vertices to a vertex in a coarser graph*/

  int 
    nedges,  ncon,
    fnedges, fnvtxs,
    cnedges, cnvtxs;
  
  int 
    *cadjncy, *cadjwgt, *cadjwgtsum, *cvwgt, *cadjidx, 
     *adjncy , *adjwgt,  *adjwgtsum,  *vwgt,  *adjidx,   
     *map, *match;  

  /*****************************/
  /* Init                      */
  /*****************************/
  fnvtxs        = graph->nvtxs;
  fnedges       = graph->nedges;
  adjwgtsum     = graph->adjwgtsum;
  adjncy        = graph->adjncy;
  adjidx        = graph->adjidx;
  adjwgt        = graph->adjwgt;
  map           = graph->map;
  match         = graph-> match;
  vwgt          = graph->vwgt;
  
  child         = LIBG_Init();

  cnvtxs        = child->nvtxs = numVtxs;
  ncon          = child->ncon  = graph->ncon;
  child->nedges = fnedges;
  
  LIBG_Allocate(child);
    
  cadjidx    = child->adjidx;
  cadjncy    = child->adjncy;
  cvwgt      = child->vwgt;
  cadjwgt    = child->adjwgt; 
  cadjwgtsum = child->adjwgtsum;
  
  cnedges = cadjidx[0] = 0; 
  
  /*Initialize the vertices list of the coarser graph (-1 => not contracted)*/
  for (i =0; i<cnvtxs; i++)
    child->vwgt[i] = -1;
  
  record  = U_mallocset(cnvtxs,int,-1,"Contract error code 2\n");
  mapadj  = U_malloc(fnedges,int,"Contract error code 3\n");

  /* map the adjacency matrix of the finer graph */
  for (i = 0; i < fnedges; i++)
    mapadj[i] = map[adjncy[i]];	  

  /*****************************/
  /* Contraction               */
  /*****************************/
  cnvtxs = 0;
  for (i = 0; i < fnvtxs; i++) {
    v = guide[i];

    /*if the label to which vertex v is mapped is equal to cnvtxs, 
      collapse it and the vertex to which it matches*/
    if (map[v] != cnvtxs)
      continue;
      
    if (ncon == 1)
      cvwgt[cnvtxs] = vwgt[v];

    /*****************************/
    /* map v in coarser graph    */
    /*****************************/
    cadjwgtsum[cnvtxs] = adjwgtsum[v];	/*Initilize adjwgsum of 
					  vertex cnvtxs in coarser 
					  graph to the adjwgtsum of 
					  vertex v in the father graph*/
    
    nedges = 0;	/*Initialize number of edges of map[v] to be zero*/
    
    /* add all edges of v in the finer graph */
    for (j = adjidx[v]; j < adjidx[v+1]; j++) {
      k = mapadj[j];	
      /* if edges is not already existing add it */
      if (record[k] == -1) {
	cadjncy[nedges] = k;
	cadjwgt[nedges] = adjwgt[j];
	record[k] = nedges++;
      /* else only update its weight */
      } else {
	cadjwgt[record[k]] += adjwgt[j];
      }
    }
    
    /**********************************************/
    /* map u (matched vertex) in coarser graph    */
    /**********************************************/
    u = match[v];	
    if (v != u) {			
      if (ncon == 1)
	cvwgt[cnvtxs] += vwgt[u];	
      
      cadjwgtsum[cnvtxs] += adjwgtsum[u];	
      
      /* add all edges of u in the finer graph */
      for (j = adjidx[u]; j < adjidx[u+1]; j++) {
	k = mapadj[j];
	/* if edges is not already existing add it */
	if (record[k] == -1) {
	  cadjncy[nedges] = k;
	  cadjwgt[nedges] = adjwgt[j];
	  record[k] = nedges++;
	/* else only update its weight */
	} else {
	  cadjwgt[record[k]] += adjwgt[j];
	}
      }
      
      /* Remove the contracted edge if existing */
      if ((j = record[cnvtxs]) != -1) {
	cadjwgtsum[cnvtxs] -= cadjwgt[j];	/* subtract its weight */
	cadjncy[j] = cadjncy[--nedges];		/* replace with last edge */
	cadjwgt[j] = cadjwgt[nedges];		
	record[cnvtxs] = -1;
      }
    }
    /*********************/
    /* clean and update  */
    /*********************/
    for (j=0; j<nedges; j++)
      record[cadjncy[j]] = -1;	
    
    cnedges += nedges;			/*add up the number of edges 
					  found in a iteration*/	
    cadjidx[++cnvtxs] = cnedges;	/*add up the edges found in a 
					  iteration to adjidx of the 
					  coarser graph*/
    
    cadjncy += nedges;		/*move the pointer preparing 
				  for receving new adjacent edges*/
    cadjwgt += nedges;		/*move the pointer preparing 
				  for receving new adjacent edges*/
    if(nedges>maxdegree)
      maxdegree = nedges;
  }
  /************/
  /* finalize */
  /************/
  child->nedges = cnedges;
  child->maxdeg = maxdegree;

  /*attach child to graph*/
  graph->child = child;
  child->father = graph;
  
  /*free tmp arrays*/
  U_free("Contract error code 3\n",&mapadj,&record,LEND);
  
  return child;
}
